corrcounts_merge <- readRDS("~/VersionControl/senescence_benchmarking/Data/corrcounts_merge.rds")
metadata_merge <- readRDS("~/VersionControl/senescence_benchmarking/Data/metadata_merge.rds")
SenescenceSignatures <- readRDS("~/VersionControl/senescence_benchmarking/CommonFiles/SenescenceSignatures_divided_newCellAge.RDS")
library(markeR)
library(ggplot2)
library(ggpubr)
library(edgeR)
Loading required package: limma
?markeR
ℹ Rendering development documentation for "markeR"
?CalculateScores
ℹ Rendering development documentation for "CalculateScores"
df_logmedian <- CalculateScores(data = corrcounts_merge, metadata = metadata_merge, method = "logmedian", gene_sets = SenescenceSignatures)
senescence_triggers_colors <- c(
"none" = "#E57373", # Soft red
"Radiation" = "#BDBDBD", # Medium gray
"DNA damage" = "#64B5F6", # Brighter blue
"Telomere shortening" = "#4FC3F7", # Vivid sky blue
"DNA demethylation" = "#BA68C8", # Rich lavender
"Oxidative stress" = "#FDD835", # Strong yellow
"Conditioned Medium" = "#F2994A", # Warm orange
"Oncogene" = "#81C784", # Medium green
"Lipid Accumulation" = "#E57373", # Coral
"Calcium influx" = "#26A69A", # Deep teal
"Plasma membrane dysruption" = "#D32F2F", # Strong salmon
"OSKM factors" = "#FFB74D", # Bright peach
"YAP KO" = "#9575CD" # Deep pastel purple
)
cellTypes_colors <- c(
"Fibroblast" = "#FF6961", # Strong Pastel Red
"Keratinocyte" = "#FFB347", # Strong Pastel Orange
"Melanocyte" = "#FFD700", # Strong Pastel Yellow
"Endothelial" = "#77DD77", # Strong Pastel Green
"Neuronal" = "#779ECB", # Strong Pastel Blue
"Mesenchymal" = "#C27BA0" # Strong Pastel Purple
)
cond_cohend <- list(A=c("Senescent"), # if no variable is defined, will be the first that appears in the ggplot
B=c("Proliferative","Quiescent"))
PlotScores(ResultsList = df_logmedian, ColorVariable = "CellType", GroupingVariable="Condition", method ="logmedian", ColorValues = cellTypes_colors, ConnectGroups=TRUE, ncol = 6, nrow = 2, widthTitle=20, y_limits = NULL, legend_nrow = 2,xlab=NULL, cond_cohend = cond_cohend)
Try scores with bidirectional signatures
bidirectsigs <- readRDS("~/VersionControl/senescence_benchmarking/CommonFiles/SenescenceSignatures_complete_newCellAge.RDS")
for (sig in names(bidirectsigs)){
sigdf <- bidirectsigs[[sig]]
sigdf <- sigdf[,1:2] # remove the third column, if applicable
if(any(sigdf[,2]=="not_reported")){
sigdf <- sigdf[,1]
bidirectsigs[[sig]] <- sigdf
next
}
sigdf[,2] <- ifelse(sigdf[,2]=="enriched",1,-1)
bidirectsigs[[sig]] <- sigdf
}
bidirectsigs
$CellAge
$CSgene
[1] "TP53" "TERF2" "MAPK14" "CDKN2A" "CDKN1A" "CCNE1" "CCNA1" "MAPKAPK5" "CBX4" "TXN"
[11] "TBX2" "STAT3" "SRF" "BMI1" "MAP2K4" "MAP2K6" "MAP2K3" "MAPK8" "MAPK3" "MAPK1"
[21] "PRKCD" "PML" "OPA1" "ATM" "MDM2" "CXCL8" "IL6" "IGFBP7" "ID1" "HRAS"
[31] "H2AFX" "POT1" "SIRT1" "KDM6B" "PLA2R1" "EZH2" "E2F3" "E2F1" "CEBPB" "CDKN2D"
[41] "CDKN2B" "CDKN1B" "CDK6" "CDK4" "CDK2" "CDC42" "RBX1" "CDC27" "CDK1" "MAML1"
[51] "CD44" "MAD2L1BP" "MAP4K4" "AIM2" "RECQL4" "ARHGAP18" "KL" "MAPKAPK2" "AURKB" "SLC16A7"
[61] "CCNE2" "HIST1H2BJ" "HIST1H3F" "CCNA2" "MCM3AP" "CDC16" "TSC22D1" "CBS" "TNFSF13" "CTNNAL1"
[71] "EED" "PNPT1" "CDC23" "RNASET2" "TP63" "CAV1" "MKNK1" "TSLP" "HIST1H2BK" "PPM1D"
[81] "HAVCR2" "CBX2" "KDM2B" "DPY30" "C2orf40" "YPEL3" "HIST2H4A" "HIST1H4L" "HIST1H4E" "HIST1H4B"
[91] "HIST1H4H" "HIST1H4C" "HIST1H4J" "HIST1H4K" "HIST1H4F" "HIST1H4D" "HIST1H4A" "HIST1H3B" "HIST1H3H" "HIST1H3J"
[101] "HIST1H3G" "HIST1H3I" "HIST1H3E" "HIST1H3C" "HIST1H3D" "HIST1H3A" "HIST2H2BE" "HIST1H2BO" "HIST1H2BC" "HIST1H2BI"
[111] "HIST1H2BH" "HIST1H2BE" "HIST1H2BF" "HIST1H2BM" "HIST1H2BN" "HIST1H2BL" "HIST1H2BG" "HIST2H2AC" "HIST2H2AA3" "HIST1H2AB"
[121] "HIST1H2AC" "HIST1H2AJ" "HIST1H4I" "HIST3H3" "CALR" "HMGA2" "PHC3" "KAT6A" "EHMT1" "SMC6"
[131] "AIMP2" "CALCA" "DEK" "MAPKAPK3" "ZNF148" "YY1" "WRN" "WNT5A" "NR1H2" "UBE3A"
[141] "UBE2E1" "UBE2D1" "UBC" "UBB" "UBA52" "CDC26P1" "TYMS" "TWIST1" "HIRA" "RPS27AP11"
[151] "HIST2H2AA4" "TP73" "TOPÂ 1,00" "TNF" "TGFB2" "TGFB1" "TFDP1" "TERT" "TERF1" "BUB1B"
[161] "BUB1" "TCF3" "TBX3" "TAGLN" "STAT6" "STAT1" "BRAF" "SREBF1" "BRCA1" "SP1"
[171] "SOX5" "SOD2" "SNAI1" "SMARCB1" "SMARCA2" "HIST2H3D" "PHC1P1" "ACD" "SKIL" "LOC649620"
[181] "SLC13A3" "LOC647654" "SMURF2" "ANAPC1" "SHC1" "CPEB1" "H3F3AP6" "ZMAT3" "RBBP4P1" "SRSF3"
[191] "SRSF1" "SATB1" "S100A6" "RXRB" "RRM2" "RRM1" "RPS27A" "RPS6KA3" "RPS6KA2" "RPS6KA1"
[201] "RPL5" "RNF2" "RIT1" "RING1" "BCL2L1" "RELA" "BCL2" "CCND1" "RBP2" "RBL2"
[211] "RBL1" "RBBP7" "RBBP4" "NTN4" "RB1" "IL21" "RAN" "RAF1" wrap_title_aux "TNRC6C"
[221] "KIAA1524" "EP400" "CNOT6" "CBX8" "PTEN" "SEPN1" "BACH1" "PSMB5" "PROX1" "PRL"
[231] "MAP2K7" "MAP2K1" "MAPK10" "MAPK9" "MAPK11" "MAPK7" "PRKDC" "RNF114" "PRKCI" "ATF7IP"
[241] "MFN1" "PRKAA2" "CDKN2AIP" "RBM38" "PRG2" "HIST2H4B" "HJURP" "TMEM140" "PBRM1" "Mar-05"
[251] "PPARG" "PPARD" "POU2F1" "TERF2IP" "ERRFI1" "H2BFS" "PLK1" "PLAUR" "PIN1" "PIM1"
[261] "PIK3CA" "PHB" "PGR" "PGD" "PIAS4" "PDGFB" "SIRT6" "ANAPC11" "ANAPC7" "ANAPC5"
[271] "WNT16" "FZR1" "ZBTB7A" "ERGIC2" "PCNA" "FIS1" "PAX3" "NOX4" "MINK1" "PEBP1"
[281] "YBX1" "NINJ1" "NFKB1" "H2AFB1" "NDN" "NCAM1" "NBN" "MYC" "MYBL2" "MSN"
[291] "ASS1" "LOC441488" "MRE11A" "MOV10" "MMP7" "MIF" "MAP3K5" "MAP3K1" "MECP2" "MCL1"
[301] "MAGEA2" "SMAD9" "SMAD7" "SMAD6" "SMAD5" "SMAD4" "SMAD3" "SMAD2" "SMAD1" "MAD2L1"
[311] "MXD1" "MIR34A" "MIR30A" "MIR299" "MIR29A" "MIR22" "MIR217" "MIR21" "MIR205" "MIR203A"
[321] "MIR191" "MIR146A" "MIR141" "MIR10B" "ARNTL" "LMNB1" "LMNA" "LGALS9" "RHOA" "KRT5"
[331] "KRAS" "KIT" "KIR2DL4" "KCNJ12" "JUN" "JAK2" "ITGB4" "IRS1" "IRF7" "IRF5"
[341] "IRF3" "ING1" "IDO1" "ILF3" "IL15" "IL12B" "CXCR2" "IL4" "IGFBP5" "IGFBP3"
[351] "IGFBP1" "IGF1R" "IGF1" "H3F3AP5" "IFNG" "IFI16" "IDH1" "ID2" "HIST2H3A" "BIRC5"
[361] "HSPB1" "HSPA9" "HSPA5" "HSPA1A" "APEX1" "HNRNPA1" "FOXA3" "FOXA2" "FOXA1" "HMGA1"
[371] "HIF1A" "ANXA5" "HELLS" "HDAC1" "H3F3B" "H3F3A" "HIST1H2BB" "HIST1H2BD" "H2AFZ" "HIST1H2AD"
[381] "HIST1H2AE" "ANAPC4" "ANAPC2" "UBN1" "SENP1" "GUCY2C" "GSK3B" "UHRF1" "BRD7" "NSMCE2"
[391] "PTRF" "GPI" "GNAO1" "RPS6KA6" "TNRC6A" "AGO2" "B3GAT1" "DNAJC2" "GJA1" "AGO1"
[401] "EHF" "TINF2" "LDLRAP1" "ULK3" "GAPDH" "ABI3BP" "ASF1A" "HIST1H2BA" "G6PD" "ACKR1"
[411] "MTOR" "CDC26" "CNOT6L" "FOS" "CABIN1" "MORC3" "SUZ12" "NPTXR" "CBX6" "SIRT3"
[421] "CRTC1" "PPP1R13B" "SUN1" "SMC5" "TNRC6B" "FOXO1" "FOXM1" "TNIK" "SCMH1" "DKKÂ 1,00"
[431] "FGFR2" "FGF2" "HEPACAM" "FANCD2" "EWSR1" "ETS2" "ETS1" "ESR2" "ERF" "AKT1"
[441] "EREG" "ERBB2" "ENG" "ELN" "CRTC2" "EIF5A" "EGR1" "EGFR" "EEF1B2" "AGO4"
[451] "AGO3" "EEF1A1" "PHC2" "PHC1" "ABCA1" "E2F2" "DUSP6" "DUSP4" "HBEGF" "AGT"
[461] "DNMT3A" "AGER" "DKC1" "DAXX" "CYP3A4" "CTSZ" "CTSD" "CSNK2A1" "E2F7" "PARP1"
[471] "HIST3H2BB" "HIST2H3C" "JDP2" "HIST4H4" "CLU" "CKB" "RASSF1" "CHEK1" "TOPBP1" "UBE2C"
[481] "KIF2C" "BTG3" "EHMT2" "GADD45G" "NEK6" "ZMYND11" "SPINT2" "CENPA" "AGR2" "CEBPG"
[491] "HYOU1" "TADA3" "MCRS1" "NDRG1" "ANAPC10" "CDKN2C" "ZMPSTE24" "PSMD14" "NAMPT" "RAD50"
[501] "TRIM10" "DNM1L" "BCL2L11"
$GOBP_CELLULAR_SENESCENCE
[1] "AKT3" "MIR543" "CDK2" "CDK6" "CDKN1A" "ZMPSTE24" "CDKN1B" "CDKN2A" "CDKN2B" "CITED2" "KAT5" "PLK2"
[13] "NEK6" "ZNF277" "CGAS" "COMP" "MAPK14" "VASH1" "PLA2R1" "SMC5" "SIRT1" "MORC3" "NUP62" "ABL1"
[25] "ULK3" "RSL1D1" "FBXO5" "FBXO4" "MAGEA2B" "NSMCE2" "H2AX" "HLA-G" "HMGA1" "HRAS" "ID2" "IGF1R"
[37] "ING2" "KIR2DL4" "ARG2" "LMNA" "BMAL1" "MIR10A" "MIR146A" "MIR17" "MIR188" "MIR217" "MIR22" "MIR34A"
[49] "MAGEA2" "MAP3K3" "MAP3K5" "MIF" "MNT" "ATM" "NPM1" "YBX1" "OPA1" "PAWR" "ABI3" "FZR1"
[61] "WNT16" "SIRT6" "PML" "PRMT6" "PRELP" "PRKCD" "MAPK8" "MAPK11" "MAPK9" "MAPK10" "MAP2K1" "MAP2K3"
[73] "MAP2K6" "MAP2K7" "B2M" "ZMIZ1" "PTEN" "MIR20B" "RBL1" "BCL6" "MAP2K4" "BMPR1A" "SPI1" "SRF"
[85] "BRCA2" "NEK4" "TBX2" "TBX3" "MIR590" "TERC" "TERF2" "TERT" "TOP2B" "TP53" "TWIST1" "WNT1"
[97] "WRN" "SMC6" "KAT6A" "ZKSCAN3" "HMGA2" "CALR" "YPEL3" "ECRG4" "MAPKAPK5" "TP63" "PNPT1" "DNAJA3"
[109] "EEF1E1" "NUAK1"
$GOBP_NEGATIVE_REGULATION_OF_CELLULAR_SENESCENCE
$GOBP_POSITIVE_REGULATION_OF_CELLULAR_SENESCENCE
$HernandezSegura
$REACTOME_CELLULAR_SENESCENCE
[1] "CDC27" "E2F2" "SCMH1" "MRE11" "MAP2K3" "MAPK9" "ANAPC4" "MAP2K4" "MAP4K4" "RPS6KA2" "UBE2D1" "EED"
[13] "MAP2K7" "TNRC6C" "MAPKAPK5" "ANAPC5" "TNRC6A" "TINF2" "AGO1" "CDC23" "CABIN1" "MAPK1" "HIRA" "TNRC6B"
[25] "E2F1" "RBBP7" "MAPK3" "ACD" "NBN" "CCNE1" "FZR1" "ERF" "CDK6" "H2AZ2" "EZH2" "MAPK8"
[37] "UBE2S" "MAP2K6" "NFKB1" "MAPK10" "ANAPC15" "CDKN1B" "PHC1" "ASF1A" "MAPK14" "E2F3" "LMNB1" "RAD50"
[49] "TFDP2" "MAPKAPK3" "IL1A" "RPS6KA1" "UBN1" "RNF2" "CDKN2C" "CDK2" "H1-3" "H1-1" "H2BC11" "CDKN1A"
[61] "ID1" "AGO3" "POT1" "CDKN2D" "CDC16" "H3-3B" "KDM6B" "TERF2" "CCNA1" "PHC2" "AGO4" "ETS1"
[73] "CDK4" "MDM2" "IL6" "TXN" "HMGA1" "RB1" "MINK1" "TP53" "ANAPC11" "CBX8" "CBX4" "RPS27A"
[85] "CCNA2" "H2BC1" "TERF1" "CDKN2B" "CDKN2A" "ATM" "HMGA2" "UBC" "VENTX" "ANAPC1" "TNIK" "MOV10"
[97] "ETS2" "H2BC5" "H4C8" "RBBP4" "MAPKAPK2" "H3-3A" "IGFBP7" "ANAPC10" "ANAPC16" "MAPK7" "TERF2IP" "H3-4"
[109] "BMI1" "H1-4" "STAT3" "CXCL8" "UBE2E1" "UBB" "FOS" "IFNB1" "CEBPB" "KAT5" "RELA" "PHC3"
[121] "CBX2" "UBE2C" "CCNE2" "ANAPC2" "CDC26" "RPS6KA3" "JUN" "SUZ12" "H2AC6" "H2BC4" "EHMT1" "EP400"
[133] "H3C13" "CBX6" "H2AC20" "H1-5" "H2BC21" "H2BC13" "MAPK11" "SP1" "H1-2" "H2AX" "H1-0" "ANAPC7"
[145] "H2AC7" "H2BC26" "H4C3" "H3C12" "H4C11" "H3C4" "MAP3K5" "H4C16" "H2BC12" "TFDP1" "MDM4" "H3C14"
[157] "H3C15" "RING1" "EHMT2" "UBA52" "H2AJ" "H4C15" "H4C14" "H4C12" "H2BC14" "H2BC8" "H3C8" "H2AB1"
[169] "H2BC6" "H4C6" "H2BC17" "H3C6" "H4C13" "H3C11" "H2BC9" "H3C1" "H4C9" "H2AC14" "H2BC3" "H4C5"
[181] "H2AC8" "H4C4" "H2BC7" "H3C7" "H2AC4" "H2BC10" "H4C1" "H4C2" "H3C10" "MIR24-2" "MIR24-1" "H3C2"
[193] "H3C3" "H2AC18" "H2AC19"
$SAUL_SEN_MAYO
$SeneQuest
NA
wrap_title_aux <- function(title, width = 30) {
if (nchar(title) <= width) {
return(title) # No need to wrap if it fits
}
wrapped_title <- ""
while (nchar(title) > width) {
# Find positions of capital letters and symbols near the wrap point
capital_pos <- gregexpr("[A-Z]", title)[[1]]
symbol_pos <- gregexpr("(_|-|:)", title)[[1]]
# Check for symbol breaks within the last few characters (width - 5 to width)
valid_symbol_breaks <- symbol_pos[symbol_pos >= (width - 5) & symbol_pos <= width]
if (length(valid_symbol_breaks) > 0) {
# If a suitable symbol is found, break at the first valid symbol
break_at <- valid_symbol_breaks[1]
} else {
# If no suitable symbol, look for capital letters within the same range
valid_capital_breaks <- capital_pos[capital_pos >= (width - 5) & capital_pos <= width]
if (length(valid_capital_breaks) > 0) {
# If a capital letter is found, break just before the capital letter
break_at <- valid_capital_breaks[1] - 1
} else {
# If no suitable symbol or capital letter, break at width
break_at <- width
}
}
# Append the wrapped line
wrapped_title <- paste0(wrapped_title, substr(title, 1, break_at), "\n")
# Update title with the remaining text after the break
title <- substr(title, break_at + 1, nchar(title))
}
# Add the remaining part of the title
wrapped_title <- paste0(wrapped_title, title)
return(wrapped_title)
}
IndividualGenes_Violins(data = corrcounts_merge, metadata = metadata_merge, genes = c("CDKN1A", "CDKN2A", "GLB1","TP53","CCL2"), GroupingVariable = "Condition", plot=T, ncol=NULL, nrow=2, divide="CellType", invert_divide=FALSE,ColorValues=senescence_triggers_colors, pointSize=2, ColorVariable="SenescentType", title="Senescence", widthTitle=16,y_limits = NULL,legend_nrow=4, xlab="Condition",colorlab="")
Using gene as id variables
CorrelationHeatmap(data=corrcounts_merge,
metadata = metadata_merge,
genes=c("CDKN1A", "CDKN2A", "GLB1","TP53","CCL2"),
separate.by = "Condition",
method = "pearson",
colorlist = list(low = "#3F4193", mid = "#F9F4AE", high = "#B44141"),
limits_colorscale = c(-1,0,1),
widthTitle = 16,
title = "test",
cluster_rows = TRUE,
cluster_columns = TRUE,
detailedresults = FALSE,
legend_position="right",
titlesize=20)
Warning: Heatmap/annotation names are duplicated: pearson's coefficient
Warning: Heatmap/annotation names are duplicated: pearson's coefficient, pearson's coefficient
Warning: `legend_height` you specified is too small, use the default minimal height.
Warning: `legend_height` you specified is too small, use the default minimal height.
Warning: `legend_height` you specified is too small, use the default minimal height.
#' @importFrom edgeR DGEList
#' @importFrom stats prcomp
#' @import ggplot2
#' @importFrom ggpubr ggarrange
plotPCA <- function(data, metadata, genes=NULL, scale=FALSE, center=TRUE, PCs=list(c(1,2)), ColorVariable=NULL,ColorValues=NULL,pointSize=5,legend_nrow=2, legend_position=c("bottom","top","right","left"),ncol=NULL, nrow=NULL){
legend_position <- match.arg(legend_position)
if (is.null(genes)){
genes <- row.names(data)
}
data <- data[row.names(data) %in% genes, , drop=F]
if (!nrow(data)>1) stop(paste0("Error: Number of genes should be >1; In your data you have only found the gene ",genes))
# Ensure metadata matches sample order if provided
if (!is.null(metadata)) {
colnames(metadata)[1] <- "Sample"
rownames(metadata) <- metadata$Sample
metadata <- metadata[colnames(data), , drop = FALSE]
y <- edgeR::DGEList( log2(data+1), samples= metadata)
} else {
y <- edgeR::DGEList( log2(data+1))
}
nPCs <- max(unlist(PCs)) # get the maximum number of PC based on the user's choice
PCAdata <- stats::prcomp(t(y$counts), scale=scale, center=center)
PCAcounts <- PCAdata$x
PCAcounts <- as.data.frame(PCAcounts)
if (nPCs > ncol(PCAcounts)) stop("Error: Number of genes too low for number of chosen PCs. Please reduce number of PCs.")
PCAcounts <- cbind(PCAcounts[,1:nPCs],y$samples)
pltList <- list()
for (pc in PCs){
pc <- unlist(pc)
ev = PCAdata$sdev^2
pc_x <- round(100*ev[pc[1]]/sum(ev),2)
pc_y <- round(100*ev[pc[2]]/sum(ev),2)
plt <- ggplot2::ggplot(PCAcounts, ggplot2::aes_string(y = paste0("PC",pc[1]), x = paste0("PC",pc[2])))
# Add jittered points, optionally colored by ColorVariable.Default: Brewer Pallette "Paired"
if (!is.null(ColorVariable)) {
plt <- plt + ggplot2::geom_point(ggplot2::aes_string(fill = ColorVariable), size = pointSize, alpha = 0.5, shape=21, color="black")
} else {
plt <- plt + ggplot2::geom_point(size = pointSize, alpha = 0.5, shape=21, color="black", fill="#D8D8D8")
}
# If ColorValues is provided, use a manual color scale; otherwise, if ColorVariable is provided,
# use a default brewer palette.
if (!is.null(ColorValues)) {
plt <- plt + ggplot2::scale_fill_manual(values = ColorValues)
} else if (!is.null(ColorVariable)) {
plt <- plt + ggplot2::scale_fill_brewer(palette = "Paired")
}
# Add axis labels (including variance)
xlab <- paste0("PC",pc[1],": ",pc_x,"% variance")
ylab <- paste0("PC",pc[2],": ",pc_y,"% variance")
titleplot <- paste0("PC",pc[1], " vs PC",pc[2])
plt <- plt + ggplot2::labs(fill = "", x = xlab, y = ylab, title=titleplot)
# Change theme
plt <- plt +
ggplot2::theme_bw()+
ggplot2::theme(axis.text.x = ggplot2::element_text(angle = 45, hjust=1),
plot.title = ggplot2::element_text(hjust = 0.5),
legend.position = "bottom")
# Adjust legend rows if legend_nrow is specified
if (!is.null(legend_nrow)) {
plt <- plt + ggplot2::guides(fill = ggplot2::guide_legend(nrow = legend_nrow, position = legend_position))
}
# Add reference lines
plt <- plt +
ggplot2::geom_vline(xintercept=0, linetype="dotted") +
ggplot2::geom_hline(yintercept=0, linetype="dotted")
pltList <- c(pltList, list(plt))
}
n <- length(pltList)
if(n==1){
plt <- pltList[[1]]
} else {
# Determine grid layout
if (is.null(ncol) && is.null(nrow)) {
ncol <- ceiling(sqrt(n))
nrow <- ceiling(n / ncol)
} else if (is.null(ncol)){
ncol <- ceiling(n / nrow)
} else if (is.null(nrow)){
nrow <- ceiling(n / ncol)
}
if (!is.null(ColorVariable)) {
plt <- ggpubr::ggarrange(plotlist = pltList, ncol = ncol, nrow = nrow, common.legend = TRUE, align = "h")
} else {
plt <- ggpubr::ggarrange(plotlist = pltList, ncol = ncol, nrow = nrow, align = "h")
}
}
print(plt)
invisible(list(plt = plt,
data = PCAcounts))
}